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Abstract 

!> ' 
1^^ ■ Using a stochastic algorithm introduced in a previous paper, we 

Q> , study the finite size volume corrections and the fluctuations of the 

ground state energy in the Sherrington-Kirkpatrick and the Edwards- 
Anderson models at zero temperature. The algorithm is based on a 
suitable annealing procedure coupled with a balanced greedy-reluctant 

00 ! strategy that drives the systems towards the deepest minimum of the 

^^ ■ energy function. 



1 Introduction 

Finding solutions of computationally hard problems is an outstanding issue 
in applied science. A classical example of hard problem is the search of the 
optimal configurations of a functional with many local minima and a com- 
plex landscape: in fact, in this case the solution is usually achieved with a 
computational effort that increases exponentially with the dimension of the 
problem. Typically, such functionals arise in the modelling of competing 
interactions among the components of large system. In some interesting 
cases in physics and in other field (for example in biology, in economy, in 
computer sciences,...), the interactions depend on some frozen-in structural 

1 



disorder [T] , and any realization of the disorder defines a sample of the sys- 
tem. Obviously, the thermodynamic properties of such systems, like spin 
glasses, fluctuate from sample to sample. Indeed, the statistical variables 
of the system (spins, in the case of spin glasses) interact via some random 
potential which models the quenched disorder. In this paper we consider 
the fluctuations, with respect to disorder realizations, of the low-lying states 
energy of spin glasses (thus the functional to be optimized here is the energy 
of the configurations of spins). Disorder- induced fluctuations are particu- 
larly relevant for the physics of systems at low temperatures. In fact, the 
properties of such systems are largely dominated by the states with min- 
imal energy into which the Gibbs-Boltzmann distribution collapses as the 
temperature is decreased. 

In order to be more specific, let us consider a disordered model with finite 
volume A^ and described by a Hamiltonian (or energy function) Hn{J, S) 
which depends on a set of statistical variables S (the N degrees of free- 
dom of the system) and on a set of variables J, representing the quenched 
disorder of the system, which are randomly sampled according to a given 
probability distribution. For a fixed disorder realization J and at zero tem- 
perature, the relevant states of the system are only the ground states, that 
is the configurations S whose energy (for degree of freedom) is given by 
e]\[ = j^ mills Hn{J, 'S')- Being the minimum of a set of random variables, 
the energy density ejsf is also a random variable whose probability density 
function (PDF) Pn{gn) follows some (a priori unknown) extreme values 
statistics: we are interested in the numerical study of this distribution in its 
large volume limit. In particular, we will mainly focus on the mean value 
^N = {s-n) and on the standard deviation aN = {{£%) — {^n)'^Y^'^ (here (. . .) 
represents the average with respect to disorder), but we will also try to 
give some insight into the shape of pAr(e7v). In order to study their scaling 
behavior, we assume that e^ and a^ have a definite limit as A^ goes to 
infinity: e^ — > e^o and ctat — > 0; in other words we suppose that extensivity 
and self- averaging of the energy density e^ hold (these properties are, in 
fact, satisfied in the models that we are going to study). More precisely, 
assuming a power-law scaling behavior for the relevant quantities: 

Cat = Coo + bN^'^ + ■■■ , aN = aN"" + ■■•,; (1) 

we aim at estimating the exponents uj and p. 

Since extensivity and self-averaging imply a trivial distribution in the large 
volume limit (pAr(x) -^ 5{x — Coo)), we inquire into the scaling behavior of 
Pn{.x) (if any) by studying the centered and scaled variable xn = {^n — 
^n)/o'n and seeking for its distribution Poo{x) in the infinite volume limit. 



The previous issues have been numerically addressed by several authors in 
recent years 0, [3] and a variety of spin glasses models have been considered 
[UlSlin]. All these studies call for large computational resources. Indeed, 
the computational effort required to solve the intrinsically complex prob- 
lem of finding the ground state of Hj<4 [1], has to be multiplied by the large 
number of realizations requested in the computation of the disorder aver- 
ages over systems with larger and larger volumes. In the cited papers [2]-[S| 
many different techniques are employed in the approximation of the ground 
state. Indeed, the relevance of solving intrinsically complex problems, has 
prompted many authors to devise efficient algorithms implementing many 
different approaches to this issue. In this study we rely on a stochastic al- 
gorithm, inspired by the classical simulated annealing technique, which was 
introduced in [7]. In that paper the algorithm was challenged against the 
search of the ground state e^ in the Sherrington-Kirkpatrick model (SK), 
which has become the standard of NP-complete problems. Our purpose 
here is to push forward the validation of the algorithm by studying in some 
details the probability distribution Pn{.^n)- Moreover, forcing the scope of 
the algorithm that was optimized in [7] for SK, we extend the study to a 
preliminary analysis of the ground state of Edwards- Anderson model (EA). 
All the results that we are going to describe are completely consistent with 
ones already present in the literature. 

The paper is organized as follows. In section 2 we introduce the spin glasses 
studied in this work, namely the Sherrington-Kirkpatrick and the Edwards- 
Anderson models, and review briefly some theoretical results that are rele- 
vant to the present issues. We sketch the algorithm in section 3, referring 
to the original paper for details. The results are given in sections 4. 

2 The models 

The Ising spin glass models we study are the Sherrington-Kirkpatrick model 
(SK) |8], whose Hamiltonian is defined by 

N 

Hn{J, S) = — 2_^ JikSiSk, (2) 

i,k-l 
i<k 

and the Edwards-Anderson model (EA) [U], whose the Hamiltonian is 

Hn{J-, S) = — 2_^ JikSiSk, (3) 

<ik> 



where Si = ±1 for i = 1 , . . . , A^ are Ising spin variables which interact 
through couphngs Jj^. For both models, Jik are independent identically 
distributed gaussian random variables, with zero mean and variance -^ for 
SK and 1 for EA. With this choice of the variance the extensivity and 
self- averaging of the ground state energy density of the two models are 
guaranteed. 

For three dimensional EA, which is a short range spin glass model, the sum 
in ([3]) is over nearest neighbor spins on a given square or cubic lattice of 
linear size L (thus the number of nodes is A^ = L'^, where d is the dimension 
of the lattice), while for SK, which is the mean field approximation of EA, 
the sum is over all the spins. (For the physical origin of these models and 
a general overview on spin glass, see [TU] ) 

The SK model in the low temperature phase was solved through the replica 
symmetry breaking ansatz (RSB) by G. Parisi [1]. This solution, which is 
universally believed to be true for SK, is still a debated issue in the math- 
ematical physics community, because rigorous proofs of some fundamental 
properties of the RSB scenario (for instance ultrametricity pQ) are lacking. 
However, recently some progress has been done in this direction with the 
mathematical proof of the Parisi formula for the free energy density [11], [12] . 
In the framework of the RSB theory, the knowledge of the ground state en- 
ergy in the thermodynamic limit (e^ = —.7633.. [1]) and the probability 
distribution of large deviations of the free energy at all temperatures [13] 
are available, but there are no exact analytic results for the limiting be- 
havior of the PDF. A computation performed at the critical temperature 
Tc and at the de Almeida-Thouless line in [T^, gives a scaling exponent 



of the internal energy density equal to |. In a more recent paper [15] a 
heuristic argument is presented which leads the same value just below T^. 
It is natural to extrapolate this result to zero temperature obtaining for the 
ground state energy density uj = ^. Indeed in [15] the authors argue, on the 
basis of the numerical evidence, that this should be the scaling exponent in 
the whole spin glass phase. Similarly, there are no analytical calculations 
for the sample-to-sample fluctuation exponent of the internal energy at zero 
temperature. However, different analytical estimates for the fluctuation ex- 
ponent of the free energy density are available: | [131 HSl US] or | [17]. As 
in the case oi u, we can expect the exponent p to be the same as that for 
the free energy density in the whole region below Tc [2] , and in particular at 
zero temperature. Unfortunately, given the closeness of the predicted val- 
ues, any numerical test at zero as well as at finite temperature can hardly 
provide a sharp distinction between two possibilities. 
For the Edwards Anderson model things are even messier because there is 



not a general consensus on the nature of its low temperature phase, and 
several scenarios (besides RSB) are conflicting [IB] . Regarding our concern, 
a relevant result was proven by Wehr and Aizenman [19]: they proved that 
for EA, as for any finite-range spin glass model in finite dimension d, the 
scaling exponent for the standard deviation a^ is p = |, namely that the 
ground-state energy variance grows linearly with the volume N = L'^. In 
fact, in [in] it is shown that the variance of any extensive quantity ^jv de- 
pending on random parameters is of the order of the volume A^. Even if this 
result seems to suggest a gaussian limiting behavior for {^ n — {^n))/VN, 
the authors of [I9] state that this, in general, should be false, though only 
very mild violations to the normal distribution are to be expected. How- 
ever, for EA, heuristic arguments as well as numerical evidence point to a 
gaussian limit for Pn{^n) [2]- The same distribution and the same estimate 
p = \ have been obtained by Aspelmeier and Moore with a replica theory 
calculation in 1201. 



3 The algorithm 

The numerical approach to the questions presented in Seel requires the 
computation of the minima of if at for a large sample of disorder realiza- 
tions, that is the search for the spin configuration which minimizes the 
energy (ground state configuration). The minimization of a certain func- 
tion depending on many discrete variables (Hamiltonian) is a combinatorial 
optimization problem. Since often many combinatorial optimization prob- 
lems are NP complete, they are tackled by constructing approximation algo- 
rithms, that run in acceptable amounts of computational time and have the 
property that final configurations are "close" to globally minimal ones (the 
metastable states approximate the ground state). The models we consider 
are often presented as the standard example of NP-problems. Indeed, the 
random sign (and strength) of the interaction generates frustration in the 
systems, i.e. the fact that in low energy configurations some of the couples 
will have unsatisfied interaction. Therefore, the global minimization can 
not be achieved simply by minimizing each local spin-to-spin interaction. 
As a consequence, the ground state of the systems is far from the standard 
ground state of ferromagnetic models, where all spins point in the same 
direction. 

Several numerical studies have tried different algorithms in the search of 
ground-state energies, for example gradient descendent [211 [221 [23] , simu- 
lated annealing [2H[25], genetic algorithms [21 [3], branch-and-bound algo- 



rithm [22] and extremal optimization [U [271 ES]- In this paper we use a 
stochastic algorithm which is the optimal one in a class introduced in [7]. 
The algorithm is required to reach the lowest possible minima as quickly 
as possible avoiding to get stuck in a local minimum which is still far from 
the deepest ones. In fact, in these complex systems disorder and frustration 
produce an energy function with a corrugated landscape, with a high multi- 
plicity of valleys (local minima) separated by high barriers (local maxima). 
Our algorithm generates a one spin-flip dynamics in the configuration space; 
indeed, starting from a randomly chosen initial condition, the algorithm ex- 
plores the space through a sequence of configurations obtained by inverting 
a single spin sign to pass from a spin configuration to successive one. The 
stochastic transition from a trajectory point to another is ruled by a prob- 
ability with an exponential density. 

Essentially, the idea is that at each step one generates a priori a random 
energy jump (variation) Aif with probability ft{AH) and then one moves 
in the direction that produces the nearest energy variation to the chosen 
jump. The algorithm stops in a local minimum that represents the best 
(sufficiently deep) encountered minimum. This algorithm is justified by 
statistical properties of the metastable states: these states are organized 
with a certain structure so that the dynamic evolution can be considered as 
the overlapping of a fast motion in the attraction basin of a local minimum 
and of a slow one with possibility of jumps between minima. The possibility 
of exceed the energy barriers, needed to escape from poor local minima, is 
obtained by introducing a sort of "external temperature" in the system, 
which enables random positive energy fiuctuations. The energy transitions 
are generated in accordance with the following PDF 

. . _ r e^i(*)^ se X < 1 1 _ 

St{x} - I ^_^^(,), gg ^ > Q , j^ + j^ - 1, 

in which the control parameter is mainly represented by the cooling rate 
A2(t) of the system and t is the time of the dynamics, i.e. the number 
of spin-fiips since the beginning of the algorithm execution. The choice of 
the exponential distribution is standard in statistical mechanics and reflects 
the Gibbs-Boltzmann equilibrium ensemble. The probability of energy in- 
creasing jumps is given by xTtY' with A2(t) increasing in t. This strategy is 
designed to model an initially hot system (A2(0) small) with high probability 
of positive moves, which is gradually quenched (the higher the temperature, 
the more likely are moves upwards and viceversa); when the system is cool 
(A2(t) large), positive fluctuations are absent and the decreasing trajecto- 
ries are forced to follow greedy-like paths (very large jumps deep into a 
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valley). With this type of algorithm we try to take advantage during the 
paths both of the greedy-like behavior and of the reluctant one (very small 
jumps and slow convergence). Obviously the performance of the algorithm 
in terms of lowest energy found and execution time, depends on the choice 
of the optimal cooling rate A2(t). In [7] the cooling rate, with an exponential 
dependence on the time A2(t) oc k~^, was optimized for the SK model. 
The large statistics of ground state values required by our analysis is ob- 
tained as follows. For a fixed disorder realization J, we performed riic in- 
dependent runs of the algorithm. Each run starts from an initial condition 
drawn at random from the uniform distribution. The (approximate) ground 
state energy eiy{J) is then obtained as the lowest energy of the metastable 
states sampled along the rijc trajectories. With the above strategy, the com- 
puted Cat (J) depends, in principle, on riic. Therefore, the reliability of the 
results can be tested, in a self consistent way, by computing the lowest en- 
ergy from larger and larger sets of initial conditions and by measuring the 
number of runs riic needed to stabilize the value within a given accuracy 
S. This stabilized value is assumed to be the ground state energy for the 
given realization of the model. In the case of SK, for instance, running the 
algorithm with the optimal cooling rate A2(t) [7], the dependence of riic on 
the volume is riic — 0.005 A^2.462 ^^^ ^ _ iq-s ^^^^ Fig{T]). The linear sys- 
tem sizes, the number of initial conditions riic and the number of disorder 
realizations rij used in our simulations are reported in Tab. [T] for both SK 
and EA model. 



4 Results 



We tested the power law ([T]) for the mean ground state energy density e^ of 
the SK model, sampling 17 different system volumes A^ between 20 and 225, 
with rij from 10^ up to 2 • 10^ disorder realizations each (see Tab. [1]). These 
values of rij have been considered appropriate because of the (relatively) 
fast convergence of the sample mean of eAr(J)'s to its limit value, as rij is 
increased. 

Fitting the data (A^, e^) to the power law Coo + bN"'^ with three free param- 
eters, we obtain the values Coo = -0.7627(±2.9 ■ 10"^), 6 = 0.7875(±7.95 ■ 
10-2), a; = 0.6820(±4.04 ■ 10"^). The influence of the subleading correc- 
tions neglected in ([1]) can be appreciated by restricting the fit to larger 
values of A^. Indeed, with A^ > 50 we have: e^o = -0.7637(±8 ■ lO-^),b = 
0.7383(±3.44 ■ 10"^), cu = 0.6617(±1.52 ■ 10"^), (see Fig. [2]). In both cases. 





N 


nj 
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1.5-10^ 


500 
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10'^ 


5-10^ 


EAsB 


3-4 


10^ 


5- 10^ 


5 


10'^ 


2-10^ 


6 


2-10^ 


10^ 



Table 1: Parameters of the simulations: linear system size, number of dis- 
order realizations nj and number of initial conditions rijc 



the estimated value of Coo is in good agreement with Parisi's analytical value 
e^ = —0.7633 [1]. The scaling exponent has been computed in previous nu- 
merical studies by Palassini [3], Bouchaud et al. [2], Katzgraber and Camp- 
bell |29]; all these results suggest u; = |. Our estimate of u with A^ > 50 
is compatible with this value, even though the uncertainty on the data are 
much larger than those presented in [2]. On the other hand, the error on the 
exponent can be lowered by fitting the data with only two free parameters 
and letting Coo = e^; in doing so, we obtain u = —0.6685 ± 2.23 ■ 10"^ (for 



A^ > 50) and to = 
u; = I we obtain b 



-0.6738 ± 9 ■ 10"=^ for A^ > 20. Moreover, fixing also 
0.7514 ± 4 ■ 10-3 for A^ > 50 and 6 = 0.7475 ± 6 ■ 10^^ 
for A^ > 20. Let us conclude this point by observing that the previous 
estimates of b are quite close to the value oi A = 0.77 ± 0.01 obtained in 
[15] fitting Coo + AN~3 to the values of the internal energy at T = 0.4. In 
[T5] the power law A^-3 has been found to describe closely the data in the 
whole spin glass phase. 

In Fig. [3]we represent in log-log scale the numerical data for a^ as a function 
of A^ together with the best numerical fit a^ = aN~P, which gives p = 
— .7597 ± 1.04 ■ 10"^. As already said in the introduction, in literature two 
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Figure 1: Scaling of the number of initial conditions per sample as a function 
of A^ in log-log scale together with the best linear fit for the SK model. 
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Figure 2: Mean ground state energy density e^v as a function of the volume 
N, together with the best three parameter numerical fit Coo + hN~'^ for the 
SK model. We obtain e^o = -.7637(±8 ■ 10"=^), b = .7383(±3.44 ■ 10"^) and 
uj = .6617(±1.52-10-^) 



values are conjectured for the scaling exponent of the sample to sample 
fluctuations: p = 3/4 and p = 5/6 [21 [31 [T71 [22]. As it is apparent from 
Figg. [3] and [H the exponent p = 3/4 (blue line) is consistent with our 
data, while p = 5/6 seems to be less likely to occur (red hne in Fig. [3] and 
red crosses in Fig. [![). This is in agreement with the results in literature 
in which all the estimates are smaller than 5/6. In particular, the value 
p = 3/4 was obtained in [3] and later confirmed in [2], but while in the 
former paper p = 5/6 is ruled out, the authors of the latter can not draw 
the same conclusion. A similar uncertainty is present also in simulations 
at finite temperature; indeed, in [15] at T = 0.4 the data for the internal 
energy are compatible with p = 5/6 but p = 3/4 is not ruled out. The 
value p = 3/4 is supported also by some recent theoretical results [131 fT5] 
suggesting that the sample-to-sample fluctuations should be proportional 
to N-^/^. 

The detection of the distribution pn{gn) for large A^ is obviously much 
more difficult than studying the two moments e^r and cr^. Such a study 
has been addressed in [21 [31 [U] via numerical simulations; here we will 
give some insight into this issue following closely the approach of [3]. The 
problem pertains to the statistics of extremes; in this theory the scaling of 
the minimum of a family of M random variables is studied for large M. 
Here the random variables are the energies H]\f{J, S) of the M = 2^ spin 
configurations S, the minimum is the ground state energy ej\f{J), and the 
scahng is studied using the variables x^ = {e^ — eN)/o'N- In the quoted 
papers the data are tested against the standard extreme values distributions 
of a family of identically distributed (i.i.d.) random variables: 1) Gumbel, 
if the individual distribution of the random variables is unbounded and 
decreases faster than any power law; 2) Fisher, if the distribution decreases 
as a power law; 3)Weibull, if the distribution has a cut-off. In the papers 
[21 [3] it is shown that none of the previous distributions describe the data. 
The same negative result was obtained [31 [15] testing the Tracey-Widom 
distribution for correlated variable. While it is not unexpected the failure of 
the Gumbel distribution in describing the data (since the energy levels that 
contribute to the ground state energy for the mean field SK model are not 
independent), it is much more surprising that the approximate behavior of 
the data can be found in the family of the generalized Gumbel distributions. 
These are the distributions the m-th smallest value in a set of i.i.d. random 
variables [30] 



gm[x) = wexp 



X — u X — u 
m 771 exp 
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Figure 3: Standard deviation ctat of the ground state energy as a function 
of A^ in log-log scale for the SK model. The lines represent the power law 
aN~P for the three values of p: 3/4, 5/6 and 0.7597 (the last one is obtained 
with a two parameter fit) 

where u, v and w are constant parameter. In his paper [3] Palassini found 
that the Gumbel distribution with m = 6 describes quite closely the ground 
state energy distribution of SK; our data support this statement. In Fig. [5] 
andOthe empirical distributions of xat for the volumes A^ = 50, 75, 100, 150, 
200 are reported together with g^lx) and the standard normal distribution. 
The values of parameters in gsix), taken from [3], are chosen imposing zero 
mean and unit variance. Even if the data for the larger sizes, see FigJHl are 
much more noisy than those for the smaller ones, the agreement with QqIx) 
seems fairly clear for all the sizes. A further support can be obtained by 
comparing the curtosis and the skewness of geix) with the same cumulants 
of the variable xn at different sizes, see table |2l 

The second family of spin system we study is the Edwards- Anderson model. 
As we said in the introduction, this should be considered as a preliminary 
study of the behavior of the algorithm in exploring the low energy configura- 
tions of the finite dimensional spin glass model. Indeed, we use the cooling 
rate optimized for the SK model, both for the two-dimensional and three di- 
mensional cases (square and cubic lattices of linear size L). However, while 
in the three dimensional case the problem of finding the ground states is 
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Figure 4: Plots of a^N^/^ and a^N^/^. 
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Figure 5: Scaling plot of the ground 
state energy PDF. The dotted line is 
the Gumbel PDF with m = 6 {u = 
.2011219, w = 2.348408, u; = 165.5589). 
All PDF {N = 50, 75, 100) are normal- 
ized to one and have zero mean and 
unit variance. 



Figure 6: Scaling plot of the ground 
state energy PDF. The dotted line is 
the Gumbel PDF with m = Q [u = 
.2011219, t; = 2.348408, m; = 165.5589). 
All PDF {N = 150, 200) are normalized 
to one and have zero mean and unit 
variance. 
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N 


Kurtosis 


Skewness 


50 


3.4249 


-.4451 


75 


3.3148 


-.4052 


100 


3.3464 


-.4046 


150 


3.5515 


-.4346 


200 


3.3312 


-.3599 


96{x) 


3.3535 


-.4247 



Table 2: Kurtosis and skewness for SK model 



NP-hard, the two dimensional one is polynomial [31]. For this reason the 
three dimensional case is only partially studied, i.e. only small linear sizes 
are considered. 

Nevertheless, our data for 3(i-model confirm the results present in the lit- 
erature [2]. All the parameters used in the numerical simulations both for 
two dimensional and three dimensional EA model are listed in Tab. [1] 
In order to validate our numerical experiment, it is useful to compare the 
numerical data for a^ with the theoretical prediction p = \- In Fig. [7] we 
show a^N^ = ai^yN] while for rf = 2 our data are compatible with size 
independence, for d = 3 a first sketchy study seems to show no analogous 
trend. In Fig. [H]we represent (for d = 2) the best numerical fit aN~P which 
gives p = 0.501 ± 6 ■ 10~^ (and a = .731 ± 1.4 ■ 10~^) in accordance with 
the expected value. This confirms the validity of our algorithm and the 
consistency of our numerical results for this model. Such a scaling low is 
expected for the central limit theorem in the cases in which the different 
terms contributing to the ground-state energy are independent. As a matter 
of fact, for EA model these terms are correlated; so this type of scaling reads 
as an indication of weak correlation. In contrast, for SK model, the variables 
that play a role in the ground state energy are sufficiently correlated to 
prevent the central limit theorem behavior. 

A further support to the hypothesis of weak correlation comes from the 
study of the limiting shape of the distribution PAr(ejv)- Because of the 
dependence of the random variables cat (J), we can not have a Gaussian 
shape for pn{,&n) at large L (see for example [19]), even though a weakly 
Gaussian-like behavior is expected (it follows from the Brout's heuristic 
argument [TU] and also from replica theory calculations [2D]). In fact, in the 
case of the two dimensional lattice, the numerical data for small volumes 
show that the PDF follows the same Gumbel distribution (with m = Q) 
already found as the limiting behavior of the SK model. However, as L is 
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Figure 7: Standard deviation of the ground state energy times the square 
root of the volume for the EA model as a function of L. 
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Figure 8: Best fit of the standard deviation of the ground state energy 
(intensive quantity) for the EA model [d = 2) as a function of L, which 
gives p = 0.501 ± 6 ■ 10^^ 
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Figure 9: Scaling plot of the ground state energy PDF for the EA model 
(d = 2). The green line is the Gumbel PDF with ?7i = 6. The red line is a 
Gaussian PDF. 

increased, we observe that the distribution moves away from the Gumbel 

and seems to approach a Gaussian limiting shape (see Fig. [9]). 

If we consider the kurtosis and the skewness of Pn{gn) we have that they 

decay with the system size (Fig. [10] and [H]). The data in Fig. [TT] suggest 

zero limiting values as A^ — > oo in accordance with the central limit theorem 

law that predicts that the skewness scales as A^~^/^ (in Fig. [11] the data 

are completely compatible with a linear convergence to zero, confirming a 

central limit-like scaling). 

Now we move to the study of the finite-size corrections to the mean energy 

density which, to the leading order, are expected to scale as N~'^. For 

EA model the analytic knowledge of eN in the thermodynamic limit is not 
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Figure 10: Kurtosis of the ground state Figure 11: Skewness of the ground 
energy distribution for EA model {d = state energy distribution for EA model 
2 and rf = 3) as a function of L. {d = 2 and (i = 3) as a function of 

available (unlike in the case of SK model). First we consider the case d = 2. 
We fit our data (3 < L < 11) to Eq. [T]where Coo, b and u are free parameters. 
We obtain (Fig. [I2]): e^c = -1.312(±2 ■ 10"^), b = 1.637(±2.94 ■ 10"^) and 
uj = 1.251(±8.5 ■ 10~^) in accordance with [2] that provides u = 1.175 ± 5 ■ 
10"'^. For d = 3 the best fit gives a mean (intensive) ground state energy 
decaying as — 1.698 + 2.077A^^^ (see Fig [T3!) . In [2] the estimated value of c^; 
is 0.967±1-10~^, even though, due to the uncertainty on the data, any value 
between 0.933 and 1 seems possible. In order to check the consistency of our 
data with this value we compute the linear fit Coo + 6A^~°'^^^, represented 
in Fig. [H which gives e^o = -1.699(±3 • lO'^), b = 1.891(±1.32 • 10"^), in 
good agreement with the values presented in [2]. 

To conclude, we presented a numerical investigations on mean field and 
finite-dimensional spin glass models that produce finite size scaling expo- 
nents consistent with the results of many works on the subject, but produced 
with different algorithms. The study of the probability distribution of the 
ground state energy showed that in the large volume limit it becomes Gaus- 
sian for the EA model, but with a finite size behavior close to a Gumbel 
distribution which, on the other hand, is the limiting distribution for the 
SK model. Further effort is necessary to produce the optimal algorithm in 
searching for the ground state of the three-dimensional EA model. 
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Figure 12: Mean ground state energy density with errorbars and the best 
fit (intensive quantity) for the EA model {d = 2) as a function of L. 
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Figure 13: Mean ground state energy Figure 14: Mean ground state energy 

density and the best fit (intensive quan- for EA model {d = 3) as a function of 

tity) for the EA model {d = 3) as a L^^-s. 
function of L. 
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